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Direct simulation of compressible 
wall-bounded turbulence 

By Gary N. Coleman 
1. Motivation and objectives 

When analyzing many turbulent flows, the effects of compressibility can be ne- 
glected. Even some relatively high-speed flows, such as boundary layers generated 
by a supersonic aircraft, produce turbulent statistics (when scaled to account for 
mean density variations) that are similar to those found for the incompressible case 
(Bradshaw 1977). There are other situations, however, in which the nonzero diver- 
gence of the turbulence leads to behavior that is fundamentally different from that 
found at constant density (Speziale Sz Sarkar 1991). Examples include flows created 
by internal combustion engines, hypersonic flight, and supersonic combustion. It is 
with instances such as these that this project is concerned. In particular, we are in- 
terested in the effects of compressibility on turbulence near a smooth solid constant- 
temperature surface; our primary objective is an increased physical understanding 
that can be used to improve turbulence models of wall-bounded compressible flows. 

With this in mind, we have begun a direct numerical simulation (DNS) study of 
turbulence in a plane channel. Because all of the relevant spatial and temporal scales 
are to be resolved, the simulations require no subgrid scale parameterization. The 
DNS code developed by Buell (1990, 1991) to study compressible plane Couette 
flow has been modified to solve the compressible Navier-Stokes equations in the 
plane channel. The channel was chosen over the Couette flow for two reasons: the 
first was to avoid using the very large streamwise domains needed to adequately 
capture the large Couette vortical structures (Lee 1990, Buell 1991); the second 
was to make use of previous experience by considering the compressible version of 
a well established case and isolate finite Mach number effects by comparing to the 
incompressible channel (Kim et al. 1987). 

The code utilizes a Fourier- Legendre spatial discretization along with a hybrid 
third-order time advance algorithm developed to optimize the range of Mach num- 
bers that may be considered (Buell 1991). The fluid is assumed to be an ideal 
gas, with constant specific heats and constant Prandtl number, and a (power-law) 
temperature dependent viscosity. Isothermal boundary conditions are used so that 
statistically stationary solutions may be obtained. The flow is driven by a body force 
rather than a mean pressure gradient to preserve streamwise homogeneity, with the 
body force f* defined so that the bulk mass flux, Q* = p*u\dy+, is constant 
(y* is the wall-normal coordinate, b+ the channel halfwidth, />* the density, and 
the streamwise velocity; “starred” quantities represent dimensional variables, and 
an overbar denotes an average over time and the streamwise and spanwise direc- 
tions). Nondimensionalized by &*, the bulk-averaged density 0* = ^ P*dy+, 
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the bulk velocity U* = and the wall temperature T the governing 

equations become: 
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The relevant nondimensional parameters are thus (a) a Mach number, M 2 = 
U 2 /'yR+T*', (b) a Reynolds number, Re = 0*£V*6*//i* ; (c) the Prandtl number, 
Pr — (d) the ratio of specific heats, 7 = C*jC*\ and (e) the viscosity 

exponent, n> where p = T n . In the above, p = p*/p * is the dynamic viscosity, 
R+ the gas constant, k+ the thermal conductivity, and the w-subscript indicates a 
value at the wall. Equations (1) - (4) are identical to those used by Buell (1991) 
except for the presence of the nondimensional forcing, = b+f* jU 2 , in ( 2 ) and 
the use of U+ in the Mach and Reynolds numbers. Another parallel to Buell (1991) 
is the separate Mach numbers in the momentum and energy equations, M and M<t 
(“d” for “dissipation”), respectively. When M = M * , a “physical” simulation is 
obtained; when they are not equal, the results can be used to differentiate between 
compressibility effects associated with the pressure gradient terms in ( 2 ) and the 
dissipation term in ( 3 ) and thus determine the relative importance of “acoustic” 
and “variable-property” influences at a given Mach number. This procedure will 
be followed in future simulations. In the next section, preliminary results from two 
M = Md runs are presented. 


2. Accomplishments 

Two cases have been simulated, one with M = Md = 1*5, the other with 
M = Md = 3; both use Re = 3000, Pr = 0.7, 7 — 1.4, and n = 0.7. The M = 1.5 
run begins by superimposing small random numbers on the laminar plane-Poiseuille 
velocity profile and uniform density and temperature fields, while the M = 3 case 
uses a mature M = 1.5 field as initial conditions. (An attempt to also begin the 
high Mach number simulation with random numbers and laminar profiles was un- 
successful since, during the initial transient, the temperature field becomes negative 
and the code stops. Lower level random disturbances could be used to avoid this 
problem but would significantly increase the time required to obtain fully developed 
turbulence — and if the disturbances are too small, there is no guarantee that a 
turbulent state will appear. To avoid a similar high -M initialization problem, Buell 
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FIGURE 1. Mass flux history, M = 1.5 (Run m01x5g). 



FIGURE 2. Wall shear histories at y = ±1: M = 1.5; M = 3; (Runs 

m01x5g and m03x0c). 


(private communication) used eigenfunctions from a linear stability analysis of the 
Couette flow to begin his simulations.) The streamwise and spanwise domain sizes 
for both cases (in units of b t ) are 4tf and 47r/3, values that are chosen to correspond 
to a DNS run for the incompressible channel. The number of collocation points 
in the streamwise x, wall-normal y, and spanwise 2 directions are respectively 70, 
90, and 40, which leads to the fully developed but marginally resolved turbulence 
desired to begin the simulations. After they have equilibrated to statistically sta- 
tionary states, the flows represented by the figures below will in the future be fully 
resolved and studied in detail. 


The M = 1.5 mass flux history is shown in Figure 1. To maintain the constant 
flux, the code adjusts the body force at each time step so that J^(&pu 1 /dt)dy & 
~p J^[du\ j dt)dy = 0 by using the Legendre quadrature routine to integrate all but 
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FIGURE 3. Mean streamwise velocity profiles: M =1.5; M = 3. Means 

in Figures 3-6 are the result of averaging over streamwise and spanwise directions 
and time; sampling period is from = 438 to 448 for Af = 1.5 and from 419 

to 429 for M = 3. 




FIGURE 4. RMS velocity fluctuation profiles for (a) M — 1-5 and (b) M = 3 runs: 

streamwise; wall-normal; and spanwise components normalized by 

(«;)„ and nZ/pl- 
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FIGURE 5. Mean density profiles: M — 1.5; M = 3. 



FIGURE 6. RMS temperature fluctuation profiles: M = 1.5; M = 3. 

the first and last terms in (2). Note that while this method is sufficiently accurate 
in practice (the value at t = 450 is within 1% of the initial flux), since we solve for 
U{ and not pu t) the constant-flux constraint cannot be enforced exactly. 

Histories of the wall shear for the M — 1.5 and M = 3 cases are contrasted in 
Figure 2. As was true for the Couette DNS (Buell 1991), an increase in M leads 
to a drop in \du\/dy\ at the walls. The higher Mach number is also responsible for 
a more pronounced flattening of the mean streamwise velocity profile as shown in 
Figure 3. Another similarity to the Couette results (cf. Figure 6 of Buell 1991) is 
revealed in Figure 4 where the near-wall maxima in the rms velocity fluctuations 
are observed to be less distinct at larger M . The mean density variation and 
temperature fluctuations are illustrated in figures 5 and 6; the large gradients of 
mean density (and therefore temperature) near the walls produce peaks in the rms 
temperature profiles, with both the density gradient and temperature fluctuations 
increasing with Mach number. 
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3. Future plans 

Having modified the code and obtained preliminary results, we now intend to 
perform a number of case studies designed to quantify the compressibility effects. 
Our immediate task is to continue the M ~ 1.5 and M = 3 runs on a fully resolved 
grid and to begin new simulations for which M and Md are not the same, as 
discussed above. Since there is a significant difference between the M = 1.5 and 
M = 3 data presented here, for the “unphysical” Md / M runs we will first choose 
M = 3 and Md = 1.5 and then M = 1.5 and Md = 3. Using the “Md ^ M” scheme 
for the Couette flow, Buell (1991) found that “pure acoustic” (nonzero turbulent 
dilatation) effects and those due to variations of mean density are both important, 
with the former influencing the wall shear and the latter the large scale structures. 
We will determine if the same is also true for the present flow and attempt to 
incorporate the findings into compressible turbulence models. 
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